submitted to ApJ 



The Effect of Poloidal Magnetic Field on Type I Planetary 
Migration: Significance of Magnetic Resonance 

Takayuki Muto 1 , Masahiro N. Machida and Shu-ichiro Inutsuka 

Department of Physics, Kyoto University, 
Kitashirakawa-oiwake-cho, Sakyo-ku, Kyoto, 606-8502, Japan 

muto@tap . scphys .kyoto-u. ac . jp 
ABSTRACT 

We study the effect of poloidal magnetic field on type I planetary migration by 
linear perturbation analysis in the shearing-sheet approximation and the analytic 
results are compared with numerical calculations. In contrast to the unmagne- 
tized case, the basic equations that describe the wake due to the planet in the disk 
allow magnetic resonances at which density perturbation diverges. In order to 
simplify the problem, we consider the case without magneto-rotational instability. 
We perform two sets of analyses: two-dimensional and three-dimensional. In two- 
dimensional analysis, we find the generalization of the torque formula previously 
known in unmagnetized case. In three-dimensional calculations, we focus on the 
disk with very strong magnetic field and derive a new analytic formula for the 
torque exerted on the planet. We find that when Alfven velocity is much larger 
than sound speed, two-dimensional torque is suppressed and three-dimensional 
modes dominate, in contrast to the unmagnetized case. 

Subject headings: MHD — planets and satellites: formation — solar system: 
formation 



1. Introduction 

Type I planetary migration is one of serious difficulties in the theory of planet formation 
and there have been a lot of work on this topic. For an unmagnetized disk, extensive linear 
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perturbation analyses have been performed since Goldreich & Tremaine's pioneering work 
(Goldreich & Tremaine, 1979). The recent result by Tanaka et al. (2002) has shown that 
the protoplanet of 5 M m located at 5AU and embedded in the minimum-mass solar nebula 



(IHayashi et al. 19851 ) will migrate inward to the central star in 8 x 10 5 years, shorter than 
the observed time scale of protoplanetary disk, 10 7 years (see e.g., Haisch et al. 2001). Their 
results are confirmed by numerical calculations (see e.g., D'Angelo et al. 2003). 

Most of analytic works on type I planetary migration have focused on unmagnetized, 
locally isothermal disks. Goldreich & Tremaine (1979) considered the two-dimensional disk, 
that is, there is no structure in the vertical direction of the disk, and derived analytic 
formulae of the torque exerted at Lindblad and corotation resonances. In order to calculate 
the migration rate, it is necessary to calculate the difference of the torque exerted on the 
disk inside and outside of the planet. This is done by Ward (1986). Later, Artymowicz 
(1993) derived, without invoking WKB approximation, a generalized formula for Lindblad 
resonances including the cutoff of the torque for high azimuthal mode number. Three- 
dimensional analysis is performed by Takeuchi & Miyama (1998), Tanaka et al. (2002), 
and Zhang & Lai (2006). Takeuchi & Miyama (1998) derived torque formulae for some 
resonances and full analytic calculation was performed by Zhang & Lai (2006). When the 
planet is embedded in the equatorial plane of a thin disk, three-dimensional modes have 
been shown to be ineffective (Tanaka et al. 2002). 

In contrast to all the above studies, magnetic fields are supposed to be present in 
protoplanetary disks. Significant mass accretion onto the central star requires an effective 
mechanism for angular momentum transfer. At present, magneto-rotational instability, or 
MRI (Balbus & Hawley 1991), is the most likely mechanism for the generation of turbulent 
viscosity in the disks. However, whether the disk is magnetically active or dead depends 
strongly on the amount of gas (and hence, dust), number density, and the properties (size 
distribution etc.) of the dust grains because the surfaces of the dust grains are very efficient 
sites of recombination (Gammie 1996, Sano et al. 2000). The detailed analysis of the 
ionization structure with the standard cosmic ray ionization rate simply predicts that the 
planet forming region (approximately between 0.1AU and 10AU from the central star) might 
be in the dead zone where gas and magnetic field do not couple. However, various effects 
can change this prediction: the growth of dust grains to larger particles decreases the total 
surface area for recombination, and hence, increases the ionization degree. Sedimentation of 
dust grains onto the disk midplane also increases the ionization degree of most of the height 
of the disk. Yet another mechanism to increase the ionization rate is also proposed, which 
possibly removes the dead zone in the standard solar nebular model (Inutsuka & Sano 2005). 



The property of planetary migration may be totally different if magnetic field is impor- 
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tant. Nelson and Papaloizou (2004) performed numerical calculations of type I migration in 
a magnetized disk and indicated that the turbulence due to MRI would result in stochastic 
torque on the planet. A simple model for the random torque was presented by Laughlin et 
al. (2004). Terquem (2003) performed the linear analysis of the torque for two-dimensional 
laminar disk with toroidal magnetic field and indicated that when stronger magnetic field 
was exerted on the disk inside the planet's orbit than the outside, inward migration might 
be halted. Fromang et al. (2005) further investigated this situation by numerical calculation 
and the numerical results showed good agreement with linear analysis. 

In this paper, we investigate the type I planetary migration when the disk is exerted 
by poloidal magnetic field, which is a complementary analysis to Terquem (2003). As a 
first step to understand the nature of migration in a magnetized disk, we perform shearing 
sheet analysis, which assumes symmetric structure inside and outside the planet's orbit, and 
calculate the torque exerted on one side of the disk. We restrict ourselves to a laminar disk, 
the case without MRI, and derive analytic formulae of torque exerted on some important 
resonances. We perform a three-dimensional calculation and our formalism is a natural 
extension of previous studies of unmagnetized cases that is developed by Goldreich and 
Tremaine (1979) and Artymowicz (1993). For two-dimensional modes, we derive an analytic 
formula which generalizes that of Artymowicz (1993). For three-dimensional modes, we 
employ WKB approximation and derive an analytic torque formula in a strong field limit. 
We show, for three-dimensional modes, that there is a divergence in perturbed density at 
certain resonances and the torque is localized at this point. We describe how to treat the 
singularity in the wave equation. We show that two-dimensional modes are suppressed by 
poloidal magnetic field and three-dimensional modes will dominate the total torque. We 
then compare the results of the linear analysis with a numerical calculation, and show good 
agreement. Type I migration in a disk with strong poloidal magnetic field may be also 
important in the formation of planets around neutron stars (e.g., Bailes et al. 1991). 

The plan of this paper is as follows. In section [21 we describe the linear analysis. In 
section [31 we describe the numerical calculation. Section H] compares the results of linear 
analysis and numerical calculation. We discuss some possible important points that are not 
covered in this analysis in section and section is for our summary. 

2. Linear Calculation of Torque 

In this section, we derive analytic formulae for the torque exerted on the disk by the 
planet by linear analysis. The backreaction exerted on the planet causes orbital migration. 
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2.1. Basic Equations 

Fo r simplicity, we consi der only a local region around the planet using shearing sheet 



model (INarayan et al. 19871 ). Although it gives the same magnitude but the opposite sign 
of the torque between the inner and outer regions of the planet so the net torque becomes 
zero, we focus on one side of the disk in this setup and simplify the problem to understand 
the effect of magnetic field. We assume that the temperature is constant and the self-gravity 
of the disk is negligible in this local region. The orbit of the protoplanet is assumed to be 
circular on the equatorial plane of the disk. We set up local Cartesian coordinates with 
origin at the protoplanet's position and the x-, y-, and z-axes are radial, azimuthal, and 
vertical direction of the disk, respectively. We use ideal MHD equations: 

^ + V-(pv) = (1) 

dv 1 1 

+ (y-V)v = —VP - W> cff - 2Q p (e z xv)- - — B x (V x B) (2) 
at p Airp 

— = Vx(vxB) (3) 

where p, v, P, tp e s, Q p , e z , and B are the gas density, velocity, gas pressure, effective 
potential including tidal force and the planet's gravitational potential, Keplerian angular 
velocity of the protoplanet, a unit vector directed to the 2-axis, and the magnetic flux 
density, respectively. We adopt an isothermal equation of state, P = c 2 p, where c is sound 
speed. The Keplerian angular velocity of the protoplanet is given by 

where G, M c , and r p are the gravitational constant, mass of the central star, and the distance 
between the protoplanet and the central star, respectively. Our calculations are normalized 
by unit time, fip 1 , unit velocity, c, and unit length, h = c/Qp. The effective potential ip e g 
in our normalization is given by, assuming a Keplerian disk, 

4>cs = ~\i 2 ~ (5) 
2 r 

where all the quantities with tilde indicate the normalized value. The first term of the right 
hand side of equation (jHJ) is composed of the gravitational potential of the central star and 
the centrifugal potential, and higher orders in x, y, and z are neglected. We also neglect the 
z-dependence of the gravitational potential of the central star for simplicity, and consider 
later the constant background density. This greatly simplifies the calculation, and we have 
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found that it does not seriously affect the results. The second term of the right hand side of 
equation (jSJ) is the gravitational potential of the protoplanet, where fn and f are normalized 
Hill radius and the distance from the center of the protoplanet respectively. The Hill radius 
is defined by r H = (Mp/3M c ) 1//3 r p , where M p is the mass of the protoplanet. 

The background disk is assumed to have no planet. The background gas flow has a 
Keplerian shear, v = — (3x/2)e y , background density is assumed to be constant, p , and 
the background magnetic field is assumed to be constant and poloidal, B = B e z . We 
denote all the background quantities with subscript zero. 

We treat the planet as a perturber on this background disk and derive the stationary 
pattern excited by the planet, as in Goldreich & Tremaine (1979). We denote perturbed 
quantities with 5, e.g., density perturbation is denoted as 5 p. We Fourier transform in £-, y-, 
and ^-directions, i.e., we shall consider the solution of the form 5p oc exp[—i(ut — k y y — k z z)]. 
Since we consider the stationary pattern, the frequency u is zero. The perturbed quantities 
are then 



where L y and L z denote the box sizes of y- and ^-directions respectively. Imposing periodic 
boundary conditions in y- and z- directions, the wave numbers in these directions are k y = 
2irn y /L y and k z = 2nn z /L z respectively, where n y and n z are integer. We shall drop the 
subscripts k y and k z of the Fourier modes unless it is ambiguous. 

We define the Lagrangian displacement £ by 




(6) 



and the inverse transformation is 





8v 



X 



ia(x)i. 



X 1 





(9) 
(10) 



where 




(11) 



Using the Lagrangian displacement, the linearized induction equations are 

SB X = ik z B £ x , 
5B y = ik z B £ y , 



(12) 
(13) 




(14) 
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The equation of continuity is 



— + ^ + ik y i y + ik z £ z = 0, (15) 
Po dx 



and the equations of motion become, using the induction equations, 

(-T 2 - 30^ + 2in p a£ y = -(c 2 + v\) ±^ + v\ (-k% - ikj^\ - * (16) 

-a 2 £, y - 2iQ p a£ x = -(c 2 + v A )ik y — - v\ (-k y k z £ z + k 2 z £ y ) - ik y ijj p , (17) 

Po 

-o 2 i z = -c 2 ik z — - ik z tp p , (18) 
Po 

where v\ = Bl/Anpo denotes the Alfven velocity of the background gas. We shall also define, 
for later convenience, the plasma j3 by c 2 /v A . 

The equations (fT5|) , (|T6|) , (|T7|) and (|T8l) are four independent equations for four variables 
5p and £. The boundary conditions to be imposed are such that wave excited propagate 
away from the planet in both inner and outer parts of the disk. 

Once the wave pattern is derived for each Fourier mode, z-component of the torque, of 
which backreaction causes the orbital migration, exerted on the disk by the planet for each 
mode is calculated by 



T k y ,k z = -2L y L z p r p k y J Im ( 



5 Pky,k z ( X ) 



Po 



^ P k y k z (a;)Gfe, (19) 



where Im denotes the imaginary part. 



2.2. Wave Propagation Property of the Disk 



We shall investigate the wave propagation property of the disk with poloidal magnetic 
field. First, we derive a wave equation from (Tl5l)- (|T8l) . The x- and ^/-components of the 
equations of motion can be written 



where f(x) is defined by 
/(*) = - 



(a 2 + 3fi 2 - v\k 2 z )i x - 2Xl v ot y = £. 
(a 2 - v\k 2 z )i y + 2iQ p aC x = ik y f, 



{(c 2 + v 2 A y-c 2 v 2 A k 2 z } 5 -P + (a 2 -v 2 A k 2 z )^ 

Po 



(20) 
(21) 



(22) 
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Equations fl20l) and (|2"T1) are the generalization of equation (10) of Goldreich & Tremaine 
(1979). Variable f(x) is related to the perturbation of total pressure 511 = c 2 Sp + B 5B z /4:7i 
by 



/(*) 



511 

Po 



p- 



(23) 



Therefore, it is a natural extension of the variable used by Goldreich & Tremaine (1979) that 
is c 2 Sp/p + ip p . Solving for £ x and £ y , 



1 



D 



(a 2 - v\k 2 z ) £ - 2Q p ak y f 
df 



-2iSl p a-±- + (a 2 + 3fi 2 _ ik j 



(24) 
(25) 



where .D is 



D 



(a 2 -^)(a 2 -^ + 3^)-4aW 



v , (26) 

This generalizes what is denoted by D in the case of unmagnetized disk, e.g., equation (12) 
of Goldreich & Tremaine (1979). Note that in the absence of magnetic field, D is a quadratic 
function of a, while this becomes a quartic function in the present situation. From equations 
(TT5"|) and (ITS]) , we finally obtain a second order ordinary differential equation which describes 
wave excitation and propagation of the disk, 

d 2 f 



dx 2 



AS + A 2 f = S, 

dx 



where 



d ^ a 2 - v A k 2 
1 dx D 



[a 2 - c 2 k 2 z )D 



+ 



2Q p ak y d 



{(c 2 + v A )a 2 — c 2 v 2 A k 2 z } (er 2 — v A k 2 ) a 2 — v 2 A 

s = "1R 

{(C 2 + V 2 A )<7 2 - C 2 v\kVl (<T 2 - V\kf) 



?7^A^^ doc 



(InD)-k 2 , 



(27) 

(28) 
(29) 

(30) 



Imposing WKB approximation, df /dx, k z f ^> k y f, this equation simplifies to 
Schrodinger type: 

<hf +V{x)f = S (31) 



dx 2 



where 



V{x) 



(a 2 - c 2 k 2 )D 



(32) 



{(c 2 + v 2 A )a 2 — c 2 v 2 A k 2 z } (cr 2 — v 2 A k 2 z ) ' 

The regions where V(x) > are wave propagation regions, and those where V(x) < 
are evanescent. The boundary between these regions, where V(x) = or V(x) = ±oo, is 
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the resonances. Figure [T] shows the appropriately normalized potential V(x) for disk with 
(3 = 0.9 and mode k y h = 0.196 and k z h = 3.14. 

There are two or three points where V(x) =0 in one side of the disk with respect to 
the planet (either x > or x < 0), depending on the value of (3. From the condition D = 0, 
we have 

.2 2 /J , 1 To2 



We call the point with positive sign Lindblad Resonance Plus (LR+) and with negative sign 
Lindblad Resonance Minus (LR— ). In the unmagnetized disk, LR+ coincides with the usual 
Lindblad resonance, a 2 = n 2 , where k is the epicycle frequency, and LR— degenerates into 
corotation point. We note that LR— exi sts only when v 2 A k 2 > 3f H, which is exactly the same 



as the stability condition against MRI (IBalbus fc Hawley 19911 ) . When LR— does not exist 



but k z 7^ 0, the corotation region becomes a wave propagation region. Since a = at the 
corotation, this indicates that there is a mode with zero frequency but non-zero wavelength, 
and therefore, in general, there is an unstable mode. Another condition for V(x) = is 

a 2 = c 2 k 2 z . (34) 

This condition does not depend on magnetic field strength. This resonance corresponds to 
that found by Takeuchi & Miyama (1998) and is named "Vertical Resonance" by Zhang & 
Lai (2006). We shall also call this point Vertical Resonance (VR) in this paper. 

There are two points in one side of the disk where V(x) diverges. One is given by 

a 2 = v 2 A k 2 . (35) 

At this point, the radial wavelength of Alfven wave becomes zero. We shall call this point 
Alfven Resonance (AR). Note that this divergence is related to what is called Alfven reso- 
nance in plasma physics (see e.g., Stix 1992). The other point where V(x) diverges is given 
by 

c 2 v 2 A k 2 



2 u A rv z (oa\ 

a = (36) 

At this point, the wavelength of the slow mode becomes zero, and this corresponds to 
"Magnetic Resonance" found by Terquem (2003) in the analysis of toroidal field. Therefore, 
we shall also call this point Magnetic Resonance (MR) in this paper. 

For two-dimensional modes, k z = 0, only LR+ exists, and the region in the vicinity of the 
corotation is evanescent region, whereas the regions further away from LR+ are propagation 
regions. The details of the wave propagation property of the disk depend on the value of 
(3 and k z , but in general, there are three propagation regions on one side of the corotation, 
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corresponding to the three wave modes of the magnetohydrodynamics. There are only two 
propagation regions in the analysis of Terquem (2003), since the mode is restricted to k z = 
and the Alfven wave with 5v z ^ is not taken into account. Figure [2] shows the wave 
propagation property of one side of the disk in the case of f3 < 1. 

With careful investigation of the properties of the resonances, it is possible to analyze 
the wave excitation, propagation, and absorption. In this paper, however, we calculate 
the torque in some restricted cases and compare the results with numerical calculation. 
We first investigate the case when k z = 0. We then calculate three dimensional modes, 
k z ^ 0, when very stro ng magnetic field is exerted. In these cases, MRI does not occur 



( IBalbus fc Hawley 199ll ). and we expect the wave pattern becomes stationary with respect 



to the planet's motion. 



2.3. Two Dimensional Mode, k z = 

We consider the two-dimensional, or k z = 0, mode. In this case it is possible to calculate 
the torque exerted on the disk without imposing WKB approximation. When k z = 0, it is 
clear from equation ffTB"]) that the fluid particles do not move along z-axis, £ 2 = 0, and the 
effect of magnetic field appears only in the pressure term. The sound speed becomes the 
phase velocity of the fast mode, c 2 + v\, and this acts as an effective sound speed. We can 
then follow the track of Artymowicz (1993) and calculate the torque exerted on the disk by 
evaluating the angular momentum flux carried by the wave at \x\ — > oo. 

The position of the effective Lindblad resonance is given by 

a 2 (x efr )-^-c 2 ^(l+r 1 ) = 0. (37) 

When 



V > 1, (38) 

Qpcky^/l + (3~ l 

the resonances inside and outside of the corotation radius are well isolated each other, and 
the torque exerted on the disk may be evaluated by the strength of the gravitational potential 
of the planet at the resonance point. The modified formula of the torque is 

T 2 d = -ir r pPoL y L z — M/ eff , (39) 

3 ft 2 + 4c 2 fc 2 (l ^2 + c2A .2 (1 + /3 -i) 

where 



M v ^ 2 + c 2 A; 2 (l + /3-i) 

#eff = -P-(X&) ~ 2k y + Vp(Zeff). (40) 
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Details of the derivation are given in Appendix [A] 

The equation (I3"9"j) generalizes the equation (56) of Artymowicz (1993) to a magnetized 
disk. Since the effective sound speed becomes faster when poloidal magnetic fields present, 
the magnitude of the torque becomes smaller. There are two kinds of cutoff mechanism of the 
torque, as Artymowicz (1993) pointed out. One is mild cutoff that comes from the coefficient 
of \l/ e ff. The other is the sharp cutoff, which is the consequence of the fact that position of 
the effective Lindblad resonance goes further away from the corotation when magnetic field 
is stronger. In the case of planetary migration, the effect of sharp cutoff is more important, 
and the torque by two-dimensional mode is strongly suppressed when (3 <^1. We emphasize 
that we have obtained the torque formula (139]) including the torque cutoff at high k y mode 
because we have not imposed WKB approximation, and hence, not neglected terms with k y . 



2.4. Three Dimensional Mode, k z ^ 0, in the Limit of Strong Magnetic Field 

In the unmagnetized disk, the contribution from three-dimensional, or k z ^ 0, mode is 



small (ITanaka et al. 20021 ). However, it is indicated in the previous section that the two- 
dimensional mode is strongly suppressed when strong magnetic field is present. Therefore 
three dimensional mode may be important in this case. We investigate three-dimensional 
mode in the limit of /3 — 0. We impose WKB approximation in this section. 

In the limit of /3 — > 0, the resonances LR+, LR-, and AR are infinitely far away from 
the corotation and we can safely neglect their contribution. The other two resonances, MR 
and VR, degenerate and the equation ([2"7]) becomes 

<hj %f = S{z), (41) 



dx 2 



where the source term S(x) is 



a 2 



S(x) = —jt—^^- (42) 



a 



c 2 k 



In this case, waves are evanescent in the vicinity of the corotation, but there is a singularity 
in the source term at the degenerate point of MR and VR, 

*mr = c 2 k 2 z . (43) 

The subscript MR denotes the quantities evaluated at this point. Note that this divergence 
originally comes from the divergence of the source term (1301) at MR. 



In order to regularize the singularity, we consider the small viscosity effective only in 
the vicinity of this point. The viscosity is effectively taken into account by adding the small 
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positive imaginary part to the frequency, 

a(x) — > a(x) + ie, 



(44) 



where e > is small positive number (see e.g., Meyer- Vernet and Sicardy 1987). Taking the 
limit of e — > 0, 



a — ck ? 



J -dx = V I ^ \ dx — m I 5d(<t — ck z )g(x)dx 



a — ck ? 



(45) 



where g(x) is an arbitrary smooth function, V denotes the principal value of the integration 
and 5d(x) is the Dirac's delta function. 



The boundary condition we impose is that the perturbation must vanish for \x\ — > oo, 



i.e. 



f(x — > oo) oc e kzX , 
f(x — > -oo) oc e kzX . 



The solution that satisfies this condition is 

1 



/(*) 



2 \Zy 



e k * u S(u)du - e k * x j e-^ u S(u)du 

+oo 



(46) 
(47) 



(48) 



Substituting S(x), the real and imaginary part of the solution are 

1 1 



f Re/ = V I" e k ^ u - x ^ p {u) 
-V j e k * {x - u) i> p {u) 

J + QO 



ck z 
1 + -^ 



2 \cr(u) — ck z a(u) + ck 2 



ck? 



1 



2 \cr(u) — ck z a(u) + ck z 



du 
du 



(49) 



and 



QVlpky 

-Imf = < 



cxp 



exp 



2ck 2 



— exp 



2ckl 

3^~2j) ^y 



2cki 



3^~2p ^y 



— e 



cxp 



2ckl 



3^~^£) ^y 



— exp 



2c^ 2 

3^~^p ^y 



x < -x M r 
-x M r <x < x MR (50) 

x > XuR 



In the limit of (3 



0, the density perturbation is given by 
bp 1 



Po 



a 2 — c 2 k 2 



a 



(51) 
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In order to calculate the torque, we need imaginary part of Sp/po which is given by 



2c/cJm 



Sp 
Po 



-it {5 D (a - ck z ) - 5 D (a + ck z )} -5- Re/ + k z ip. 



a- 



a 2 

+ — 

vi 



V 



V 



a — ck z a + ck z 



Ivaf. 



(52) 



It is possible to show that / does not diverge at the resonance (see Appendix [B] for details) , 
and we can neglect the term with o 2 jv\ in right hand side of the equation fl5T|) when magnetic 
field is strong enough. Quantitatively, we can neglect these terms when 



cr 2 ck 2 

^ A ^pky 



(53) 



since / ~ 0(ck 2 ifj p /^l p k y ). Since a ~ c k 2 in the vicinity of the resonance, we obtain 

O h 

ckj 

When this condition is satisfied, the imaginary part of the density perturbation is 



(54) 



Im 



Sp 
LPo. 



tk z 
' 2c 



{5 D (a - ck z ) - 5 D (a + ck z )}ip p 



(55) 



The first delta function indicates the torque exerted on the outer disk, and the second inner 
disk. These torques are of the same magnitude but different in sign. The magnitude of the 
torque on one side of the disk is then, from (Tl9l) . 



T - 2W T I 

J- MR — ~3~ 



>i' J z c '• />.MH- 



(56) 



3. Numerical Calculation 

We have performed numerical calculations in order to investigate how well the equations 
( |39l) and ( |56i) describe the realistic value of the torque. We have done two sets of runs. One 
is for a two-dimensional disk. The other is for a three-dimensional thick disk. 



3.1. Numerical Methods 



We adopt the nested grid method (see, e.g., Machida et al. 2005, Matsumoto & Hanawa, 
2003) to obtain high spatial resolution near the planet. Each level of rectangular grid has 
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the same number of cells (= 64 x 256) for 2D run, while (= 64 x 256 x 16) for 3D run. 
The cell width As (I) depends on the grid level I. The cell width is reduced by half with 
increasing grid level (I — > I + 1). We use 4 grid levels (I =1,2 ■ • ■ 4) for 2D run and 5 levels 
for 3D run. The box size of the coarsest grid I = 1 is chosen to (L x ,L y ) = (64/i, 256/i) 
for 2D run and (L x , L y , L z /2) = (64/i, 256/i, 16h) for 3D run. Note that in ^-direction, the 
simulation box extends from midplane to z = L z /2. The box size of the finest grid is 
(x, y) = (2h, 8h) for 2D run and (x, y, z) = (2h, 8h, h) for 3D run. The cell width of the 
coarsest grid is As(l) = h, while that of the finest grid has As(4) = 0.125/i for 2D run and 
As(5) = 0.0625/i for 3D run. We assume the fixed boundary condition in the x-direction 
and periodic boundary condition in the ^/-direction. For z-direction, we impose a periodic 
boundary condition between z = —L z /2 and z = L z /2. 

For two dimensional calculation, we neglect the z-dependence of the planet potential, 
i.e. we adopt the potential of the form 

v x + y 

We employ the softening in the gravitational potential as follows. The gravitational 
force F by the planet is given by 

GM P 
[r + e) 

where e is the softening length, r is the distance from the planet's position, and x is the 
position vector. We choose e such that this equals the mesh size of the finest grid, i.e., 
e = 0.125/i for 2D run and e = 0.0625/i for 3D run. 

We fix the planet mass to be fn = 0.3, corresponding to 3M e planets when M c = M & 
and h/r p = 0.05, and vary the initial strength of the poloidal magnetic field. We performed 
the calculations for (3 = 00, 100, 10, 2, 0.3, 0.1, 0.01, and 0.001. 



F=7^V. (58) 



4. Comparison between Numerical Calculation and Linear Analysis 

For all two-dimensional calculations and for three-dimensional calculation with (3 = 00, 
0.01, and 0.001, we do not observe MRI and steady states are realized. This is consistent 
with the stability criterion of MRI derived from linear analysis. 

We then Fourier transform the density pattern of the steady state in y- and z-directions 
and calculate the torque exerted on one side of the disk by equation (Tl9|) . and this torque is 
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compared with the results of linear analysis. The normalization of the torque is taken to be 

T 

T = T^HHP, (59) 
We make use of FFT (e.g., Press et al. 1992). The wavenumber k we evaluate is given 

by 

27m . . 

k = —j—, (60) 

where L is the box size of the y- or z-directions and n is an integer with —N/2 < n < N/2 
where N is the mesh number. We also Fourier transform the gravitational potential of the 
planet numerically to obtain the value of ip p . 



4.1. Two-dimensional Calculation 

Figure [3] shows the stationary pattern of density perturbation obtained by two dimen- 
sional calculations for (3 = 0.01, 2, and 100. It is clear that, with increasing magnetic field, 
the amplitude of the wave becomes small and the point where waves are excited goes further 
away from the planet. Figure H] shows the torque calculated as a result of numerical calcu- 
lation for various magnetic field strength, or different j3. It is clear that the torque becomes 
weaker as the magnetic field is stronger. 

We show in figure the comparison between the results of numerical calculation and 
linear analysis, equation fl39l) . It is clear that for modes that satisfy condition fl38l) . which we 
expect that equation fl39l) gives a good approximation for the torque, numerical calculation 
and linear analysis indeed show reasonably good agreement, at least an order of magnitude, 
even though equation (159~j) estimates the torque by the value of density perturbation only at 
the position of effective Lindblad resonance. Therefore, equation (1391) is useful for estimating 
two-dimensional torque when poloidal magnetic field is exerted on the disk. 

We also checked that the numerical calculation and linear analysis are in good agreement 
for other values of (3 except for (3 = 0.001. For (3 = 0.001, since the amplitude of density 
perturbation is very small, numerical torque is dominated by small noise in the disk. 



4.2. Three-dimensional Calculation 



For (3 = 0.01 and (3 = 0.001 models of the three-dimensional calculations, we do not 
observe MRI and steady state is realized. For other parameters, we observe the instability. 
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Since we investigate the stationary pattern, we focus on results in which we do not observe 
MRI. We show in figure [6] the torque that is derived from numerical calculation for n z = 
0,1, and 2, where n z is the mode number of z-direction. It is clear that n z — 1 modes 
overwhelm the two-dimensional modes in these models. 

Figure [7] compares the torque calculated from the three-dimensional numerical calcula- 
tions and that calculated from linear analysis of n z = 1 modes, torque formula (I56I) . From 
the derivation of formula (1561) . this expression of the torque is valid when WKB condition 
k y <C k z and strong magnetic field condition (I54p are both satisfied. In the present parame- 
ter, the WKB condition is more restrictive. Since k z h = 27r/32 = 0.196 for n z — 1 mode, we 
expect that for k y h greater than this value, equation fl56|) does not give a good approxima- 
tion for the torque. Nevertheless, the result of the numerical calculation indicates that the 
equation fl5T)l) shows a very good agreement even in the modes with k y h greater than this 
limit. 

We also find that the imaginary part of the Fourier components of density perturbation 
diverges around MR, as expected from linear analysis. Figure [H] shows the profile of the imag- 
inary part of the density perturbation of (3 = 0.001 calculation for (k y h, k z h) = (0.498, 0.785). 
The position of magnetic resonance is indicated by an arrow. It is clear that density pertur- 
bation diverges at the resonance position and the contribution of the torque mostly comes 
from this divergence. The torque is localized at the magnetic resonances since waves can- 
not propagate on the disk, and, therefore, the analytic torque formula ([55]) gives a good 
approximation of the total torque, even if we consider the regions only in the vicinity of the 
resonance. 



5. Discussion 

5.1. The Strength of Three-dimensional Modes in a Thin Disk 

Tanaka et al. (2002) has shown that in the unmagnetized disk, three-dimensional modes 
are subdominant. In contrast, when poloidal magnetic field is exerted on the disk, it is 
indicated that three-dimensional modes can dominate the torque when the magnetic field is 
sufficiently strong. In this section, we briefly discuss the critical value of (3 at which k z ^ 
modes dominate the total torque in a thin magnetized disk according to the results of linear 
analysis. By "thin disk" , we refer to the disk with small aspect ratio, smaller than that we 
have used in the numerical calculation, but not two-dimensional. 

We calculate the torque for n z = modes by fl39l) and n z — 1 modes by fl56|) . The 
Fourier transformation of the planet's gravitational potential is done numerically with the 
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box size L x = 32h, L y = 128/i, L z /2 = 2h and the spatial resolution 0.125/i. We have checked 
that n z = 2 modes are smaller than n z — 1 modes. Figure [9] shows the torque for (3 = 0.01, 2, 
and 100. The vertical wavenumber of n z — 1 modes is k z h = 1.57, and WKB approximation 
is valid for k y h less than this value. For thin disk case, it is indicated, just as thick disk case, 
that two-dimensional modes are dominant for a weak magnetic field, while three-dimensional 
modes are important for a strong magnetic field case. We have investigated other values of 
(3 and it is indicated that for the disk with (3 <^ 0.1, three-dimensional modes are more 
important than two-dimensional modes. 

Since three-dimensional torque formula fl56|) is valid only in the strong magnetic field 
limit, it is not possible to extrapolate this to the case with j3 ~ 1. However, since the torque 
formula for two-dimensional modes f[59"j) does not have any restriction, we can safely conclude 
that k z = modes are always suppressed for strong magnetic field. Therefore, qualitatively, 
we expect that two-dimensional modes are suppressed for (3 ^ 1. To verify this conjecture 
quantitatively, we need careful analysis for (3^0 case, which will be presented elsewhere 
(T. Muto and S. Inutsuka 2008, in preparation). 



5.2. The Relation between Magnetic Field Strength and the Differential 

Torque 

For an unmagnetized disk, it is known that the outer torque that is exerted by the 
disk outside the planet wins over the inner torque exerte d by the di sk inside, when the 



disk gas density is larger in the inner disk than the outer (IWard 19971 ). This is the result 
of the competition of two effects. On one hand, since the inner density is larger than the 
outer, the inner torque becomes larger than the outer. On the other hand, the effect called 
pressure buffer enhances the outer torque. Considering the background disk structure, gas 
is slightly sub-Keplerian resulting from the outward pressure gradient. Since the planet is in 
Keplerian rotation, the corotation point locates slightly inside the planet, which makes the 
outer Lindblad resonance slightly closer to the planet. Calculating the difference of these 
two competing effects, the outer torque is larger than the inner torque. 

Let us now qualitatively discuss the differential torque in the disk with poloidal magnetic 
field. First, we consider the two-dimensional mode when magnetic field exerted on the disk 
inside the planet's orbit is larger than the outside. In the case without variation in density 
and temperature, the mild cutoff of the torque by magnetic field makes the outer torque 
stronger. If the disk has radially decreasing magnetic pressure distribution, the planet locates 
slightly outside the corotation point since outward magnetic pressure is exerted on the gas. 
This also enhances the outer torque, since the outer effective Lindblad resonance locates 
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closer to the planet. The outer torque is, therefore, expected to be stronger because of these 
two effects, in contrast to the toroidal magnetic field case of Terquem (2003) where the inner 
torque is stronger. This might indicate that the differential torque may be very sensitive to 
the configuration of magnetic field near the planet. 

We now turn to the magnetic resonances of k z ^ modes, effective for low (3. Since the 
formula (1561) is for the limiting case of (3 — > 0, the torque does not depend on the magnetic 
field strength. We shall propose a simple torque formula for MR which generalizes equation 
(I56p and discuss the effect of magnetic field. Firstly, we note that the relation between / 
and Sp/po is given by 

'i-wrvh^y'-^-*™- (61) 

Since / in the right hand side can be neglected when (3 = 0, we expect this term can be 
neglected even in (3 ^ 0, provided that (3 is sufficiently small. This equation also indicates 
that there is a 5-function-like divergence at MR. Neglecting the term with /, we obtain the 
following torque formula at MR for low f3, 



_ 2tt p r p k z 2 



T T ^° P z /2 (co\ 



where the value of gravitational potential is evaluated at MR. The smaller the magnetic 
field strength is, the closer towards the planet the MR position locates, which makes the 
gravitational potential at MR, Vp,mr> stronger. However, the coefficient, (1 + /3)~ 3//2 , be- 
comes smaller, which makes the evaluation complicated. Using the parameters with three- 
dimensional torque calculation, we evaluate the Fourier transform of the gravitational po- 
tential and calculate the torque. Figure [10] shows the torque calculated from the modified 
formula (1621) for f3 = 0.001, 0.01, and 0.1. It is indicated that MR torque becomes smaller 
for weaker magnetic field strength. Therefore, we expect that when the inner magnetic field 
is stronger than the outer magnetic field and the field strength is high enough for k z ^ 
modes to be dominant, the inner torque wins over outer torque, in analogous to the results 
of the analysis of toroidal field by Terquem (2003). 

When magnetic field is very strong, equation (R)2"|) indicates that the value of the torque 
is not sensitive to the strength of the field, and the differential torque can be very small. 
We consider, then, the effect of the gradient of sound speed, which changes the location of 
the resonance even in (3 — > limit. Let us consider the disk with higher sound speed inside. 
The outer MR, which nearly degenerates with VR, is closer to the planet, giving a larger 
value of gravitational potential. The coefficient of the gravitational potential in equation 
( I56p is inversely proportional to the sound speed because L z k z is the mode number of the 
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z-direction that is indifferent to the value of c. Therefore, the coefficient is smaller for the 
inner MR than the outer MR. The outer torque is expected to be stronger than the inner 
torque when magnetic field is very strong and when there is a negative gradient of the sound 
speed. 

The qualitative dependence on magnetic field of the differential torque for laminar modes 
may be summarized as follows. Consider the case where inner magnetic field is stronger than 
the outer. When magnetic field strength is weak and two-dimensional modes are dominant, 
the outer torque is more enhanced and the migration is inward. When magnetic field is 
strong enough for three-dimensional modes to be dominant, the migration can be outward. 
Note, however, that the rate and directions of migration may depend sensitively on the 
gradient of the sound speed c. Negative gradient of the sound speed may cause the inward 
migration for very low (3. 

It seems difficult to halt the inward migration in a disk with strong poloidal magnetic 
field, since outward migration may require positive gradient of sound speed. Note, however, 
that the typical magnitude of one-sided torque is always smaller than the unmagnetized 
case, as shown in figure M In the disk with (3 = 100, the torque is dominated by two- 
dimensional modes and its magnitude is approximately 10 -3 in our normalization, while in 
(3 = 0.01 case, the magnitude is smaller by about two orders of magnitude. Therefore, the 
strong magnetization of disk is expected to slow down the migration. Actually this outcome 
is analogous to the effect of increasing gas temperature, and hence, thermal pressure and 
sound speed in the disk without magnetic field. 



5.3. Comparison with Toroidal Field Case 

In this paper, we have considered a protoplanetary disk threaded by a poloidal magnetic 
field. In the analysis of the torque at the magnetic resonance, we have considered a strong 
magnetic field case and derived a torque formula (|56|) . We now consider briefly the more 
general case when toroidal component of magnetic field also exists. Although it is necessary 
to make a full, rigorous calculation including both poloidal and toroidal components of 
magnetic field in the background disk, we make a qualitative discussion by comparing the 
effect of magnetic resonances of purely toroidal and that of purely poloidal case. 

In the case of toroidal magnetic field, there is a magnetic resonance in two-dimensional 
modes too (Terquem 2003). In the strong field limit, the position of the magnetic resonance 
is given by 

2 

£MR,toroidal = ^H. (63) 



-19- 



Since the magnetic resonance in poloidal case is located at 

2 k 

^MR,poloidal = 7 7"^, (64) 
O rCy 

and in a thin disk, modes k z ^> k y are important, the magnetic resonance of a toroidal field is 
closer to the planet than the toroidal case. We also note that the Fourier components of three- 
dimensional modes of gravitational potential are much smaller than the two-dimensional 
modes, as shown by Tanaka et al (2002). Therefore, in a standard case of a planet on a 
circular orbit embedded in a disk midplane, we expect that when the net magnetic field is 
dominated by toroidal components, the effect of poloidal magnetic field is small compared 
to that of a toroidal field. A possible exception is provided by a large inclination of planet's 
orbit. In this case, the gravitational potential of planet has a large z-component (three- 
dimensional modes), and thus, the poloidal field would be important. Note, however, that 
we need a new set of analyses for the case of a planet with inclined orbit (for unmagnetized 
disks, see Tanaka & Ward 2004). 



6. Summary 

We have performed linear perturbation analysis to calculate the torque exerted on the 
planet embedded in the non-turbulent disk with poloidal magnetic field using the shearing 
sheet approximation. We have derived a second order ordinary differential equation de- 
scribing the excitation and propagation of the wave, equation (l2"7|) . and derived analytic 
expressions of the torque for two limiting cases. Equation (139]) gives the torque for two- 
dimensional modes. Equation gives the torque for three-dimensional modes for < 1 
under WKB approximation. We have compared the result of the linear analysis and numer- 
ical calculation, and found that both formulae show reasonable agreement. We have shown 
that the two-dimensional modes are suppressed when magnetic field is strong, in contrast to 
Terquem (2003) analysis of toroidal field, indicating that the property of planetary migra- 
tion may be sensitive to the configuration of the magnetic field around the planet. It is also 
indicated that when magnetic field is very strong, three-dimensional modes are more effec- 
tive than the two-dimensional modes, in contrast to the analysis of the three-dimensional 
calculation of unmagnetized disk by Tanaka et al. (2002) 

Since we have been using the shearing sheet approximation and have derived the torque 
formulae only in restricted cases, the analysis of more general cases and other resonances 
is necessary, which will appear elsewhere (T. Muto and S. Inutsuka 2008, in preparation). 
Although the equation fl56l) agrees well with the numerical calculations, this is derived under 
WKB approximation. Therefore, we need careful analysis for high k y for three-dimensional 
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modes. We also need more quantitative analysis of the differential torque, which determines 
the direction and rate of the migration of the planet. 
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A. Derivation of Torque Formula for k z = Modes 



In this section, we derive the equation (I39p . following Artymowicz (1993) formalism. 
The linearized equations of continuity and motion for k z = modes are 



— lO h — 0V X + tkyOVy = 0, 

po dx 

-ia5v x + c 2 (l + P' 1 )-^-— - 2Q p 5v y = — ^ip p , 

ax po ax 

—iaSvy + -Vt p 5v x + ik y c 2 (l + fi^ 1 )— = —ik y ip p . 
2 Po 



From these, we obtain the equation for vorticity: 

dx 2 po 



0. 



(Al) 
(A2) 
(A3) 

(A4) 



Using the equations of motion, we finally obtain the Schrodinger-type second-order ordinary 
differential equation for 5v y : 



dx 2 y+ c 2 (l + ^- 1 ) 



— y-tf p -c 2 (i + p-i)kl] SVy 



2 p llx~ y ^ p 



c 2 (l + (3- r 

The position of the effective Lindblad resonance is given by 

a 2 = n 2 p + c 2 (i + p- 1 K 



(A5) 



(A6) 
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Waves are evanescent in the region in the vicinity of the planet and the regions further away 
from the resonance are propagation regions. Waves are excited at the effective Lindblad 
resonances and propagate away from the planet to \x\ — > oo. 

We shall calculate the angular momentum flux at infinity. The angular momentum flux 
is calculated by 

F A = 2r p p L y L z Re [5vjv*] . (A7) 
Since the gravitational potential of the planet vanishes at infinity, the flux at the infinity is 



r r 4C 2 (1 + (3- l )k v 

F A {x - oo) = 2r pPo L y L z hn 



5v *Tx 5v » 



(Ah 



The solution of the wave equation (1A5I) is given by parabolic cylinder functions. Here, for 
simplicity, we assume the two resonances, inside and outside the planet, are isolated each 
other. The equation in the vicinity of the resonance is then given by 



where 



d 2 

—5^ + 27(2 - i)8v y = -S, 



(A9) 



3 



X 



7 



, _ 2 fl 2 p + c 2 (l+(3- 1 )kl 
"3 Qpcky^l + p- 1 



s 



(hpp 
dz 



l8ck y 



i> 9 



off 



The subscript "eff" denotes the quantity evaluated at the resonance. 



(A10) 
(All) 
(A12) 



We impose the boundary condition as follows. In the evanescent region, the solution 
does not grow exponentially and in the propagation region, the waves propagate away from 
the planet. The solution is then 



5v,„ 



7lS 



(2 7 )V3 



{Gi [-(2 7 ) 1 / 3 (^ - 7 )] + *Ai [-(2 7 ) 1/3 (^ - 7)] } 



(A13) 



Where Ai represents the Airy function and Gi is the solution of the equation 
( lAbramowitz fc Stegun 19701 ) 



dx 2 



Gi(x) — xGi(x) 



1 

7T 



(A14) 
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Substituting the equation (1A13I) into (IA8I) . we obtain the torque formula fl39l) . 

The condition of the isolation of the resonances is satisfied when the solution at one 
resonance does not affect the other resonance. The distance between the resonances is 
6z ies ~ 7, while the scale that the solution in the vicinity of one resonance changes in the 
evanescent region is 5z wave ~ 7 -1 / 3 . Hence, resonances are well separated each other when 

5z Tes > fc w ave, Or 7 > 1. 



B. The Evaluation of the Magnitude of / 



In this section, we evaluate the integral in the equation (T49I) . For simplicity, we set the 
planet's gravitational potential ip p to be constant. The part of the equation (1491) . 



d ue kz(u-x) 



(Bl) 



is finite. We consider the rest. Since we are working on the local Cartesian coordinate where 
inside and outside the planet are symmetric, we assume x > without loss of generality. Let 
/ be 



I = V 



-V 



— oo 

X 



due k ^ u ~ x) 
due kAx -' a) 



a(u) — ck z a(u) + ck- 
1 1 



a(u) — ck z a(u) + ck z 



then 



/ = I A - I B - I c + I D , 



I a = exp 
I B = exp 
I c = exp 
I D = exp 



-k z x 



k z ( x + 



3^~^.p ^^y 

2ck z 



Ei ( k z x 



3^2p ky 

Ei k z x + 



3i~^^ ^^y 



x — 



k z \ x + 



2ck z 



Ei -k z x + 



3 p ky y 



Ei ( — k 7 x 



3 ^^ip \\jy 

2ck 2 z 
3 ^^iip ky 



where Ei denotes the exponential integral, 

Ei(w) = V 



dt—. 
t 



(B2) 

(B3) 
(B4) 

(B5) 

(B6) 

(B7) 
(B8) 

(B9) 
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We can check I a, Ib, Ic, and Id do not diverge at infinity by the asymptotic expansion 
of the exponential integral, 

rw t i 

-dt (Bf 0) 



<x> 



t W 



and the inequality derived from the series expansion of the exponential integral, 

00 n 

EH) 
— - < 7 + lnw + e™. (Bll) 
nn\ 

n=l 

For x -> 00, I ~ 0(1). 

We now consider the vicinity of the resonance, x ~ 2ck z /3Q p k y . Although I a and Ic 
are divergent logarithmically, the combination I a — Ic does not diverge: 

pk z 5 t 

I A -Ic~V I dt- (B12) 

J-k z S 1 

where we set x ~ 2ck z /3Q p k y + 5. It is easy to show Ib,Id ~ an d therefore, / is 

order of unity for all x. Although ip p is not strictly a constant, it is smooth in the vicinity 
of the resonance. Therefore, we can safely assume this to be constant when we discuss the 
divergence at the resonance. 

In summary, from (1491) , the order of magnitude of Re/ is 

y i bp rCy J 

It is clear from the equation ([50]), the imaginary part of / is also of the same order. 
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potential for k y h=0.196, k z h=3.14, p=0.9 
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Fig. 1. — Potential V(x) given by equation ( |32i) for /3 = 0.9, fc^/i = 0.196, and k z h = 3.14. 
Resonance positions in the outer disk (x > 0) are indicated. LR+ and LR— denote Lindblad 
resonances, AR denotes Alfven resonance, VR denotes vertical resonance, and MR denotes 
magnetic resonance. The grey regions correspond to the evanescent regions. Note that 
regions \x/h\ > 15 are all propagation regions. 
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Fig. 2. — Wave propagation property for j3 < 1 disk. The parameters used in figure [Q 
corresponds to case (c). 
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Fig. 3. — Density profile obtained by the two-dimensional numerical calculation. The pattern 
of density perturbation 5p/p is indicated by false color. It is clear that the stronger the 
magnetic field, the further the point where waves are excited and the smaller the amplitude. 
Note that color scales are different for different values of (3. The x- and y-axes correspond 
to the axes of shearing-sheet, normalized by the scale height c/Q p . The elapsed time t is 
normalized by the planet's Kepler time O" 1 . Four different levels of nested grid are super- 
imposed. 
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Torque obtained by 2D numerical simulation 




k y h 

Fig. 4. — The dependence of two-dimensional torque (k z = modes) on the strength of 
magnetic field obtained by numerical calculations. The models with f3 = oo (no magnetic 
field, plus), 10 (cross), 2 (open square), and 0.1 (filled square) are shown. The torque is cut 
off drastically for models with (3 < 1. The horizontal axis denotes k y h and the vertical axis 
denotes normalized torque. 
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Fig. 5. — Comparison of the torque obtained by the two-dimensional numerical calculation 
(plus) and the linear analysis (line), the equation (|39|) for (3 = 100 (top left), (3 = 2 (top 
right), and /? = 0.1 (bottom). The horizontal axis denotes the azimuthal mode number and 
the vertical axis denotes normalized torque. 



-31 - 




Fig. 6. — /^-dependence of the torque obtained by three-dimensional numerical calculations. 
The left panel shows f3 = 0.01 and the right /3 = 0.001. The horizontal axis denotes the 
azimuthal mode number and the vertical axis denotes normalized torque. Two-dimensional 
modes are denoted by plus, three-dimensional modes with n z = 1 by cross, and n z = 2 by 
filled square. Three-dimensional modes with n z — 1 dominate the two-dimensional modes 
when the magnetic field is sufficiently strong (f3 <^ 1). 
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Fig. 7. — The comparison of n z = 1 mode torque between the analytic formula (line), 
equation (1391) . and three-dimensional numerical calculation for /3 = 0.001 (plus) and j3 = 0.01 
(cross). The horizontal axis shows the azimuthal mode number and the vertical axis shows 
the normalized torque. 
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k y h = 0.498 




5 10 15 20 



Fig. 8. — The radial profile of lm(5p/po) for the Fourier component (k y h,k z h) = 
(0.498,0.785) obtained by f3 = 0.001 numerical calculation. The profile of the torque on 
the disk depends on the imaginary part of the density perturbation [see equation ffl9l ]. 
The horizontal axis shows the radial coordinate and the arrow indicate the position of the 
magnetic resonances calculated by the linear analysis. 
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Fig. 9. — Comparison between two-dimensional torque formula (1391) and three-dimensional 
formula fl56l) for thin disk with (3 = 0.01 (top left), (3 = 2 (top right), and (3 = 100 (bottom). 
Two-dimensional torque is denoted by the solid line and three-dimensional torque is denoted 
by the dashed line. The horizontal axis denotes azimuthal mode number and the vertical 
axis denotes normalized torque. 
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Torque at Magnetic Resonance, n z =1 
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Fig. 10. — Torque at the magnetic resonance obtained by equation ( 1621) . Magnetic field 
strength corresponds to (3 — 0.001 (solid line), (3 = 0.01 (dashed line), and (3 = 0.1 (dotted 
line). The horizontal axis denotes azimuthal mode number and the vertical axis denotes 
normalized torque. 



